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Abstract 


The effect of micro-gravity on the musculoskeletal system has been well studied. Sig- 
nificant changes in bone and muscle have been shown after long term space flight. Similar 
changes have been demonstrated due to bed rest. Bone demineralization is particularly pro- 
found in weight bearing bones. Much of the current techniques to monitor bone condition use 
bone mass measurements. However bone mass measurements are not reliable to distinguish 
Osteoporotic and Normal subjects. It has been shown that the overlap between normals and 
osteoporosis is found for all of the bone mass measurement technologies: single and dual 
photon absorptiometry, quantitative computed tomography and direct measurement of bone 
area/ volume on biopsy as well as radiogrammetry. A similar discordance is noted in the fact 
that it has not been regularly possible to find the expected correlation between severity of 
osteoporosis and degree of bone loss. 

Structural parameters such as trabecular connectivity have been proposed as features 
for assessing bone conditions. In this report, we use fractal analysis to characterize bone 
structure. We show that the fractal dimension computed with MRI images and X-Ray 
images of the patella are the same. Preliminary experimental results show that the fractal 
dimension computed from MRI images of vertebra of human subjects before bedrest is higher 
than during bedrest. 
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1 Introduction 


The effect of micro-gravity on the musculoskeletal system has been well studied. Significant 
changes in bone and muscle have been shown after long term space flight. Similar changes 
have been demonstrated due to bed rest. Bone demineralization is particularly profound 
in weight bearing bones. Much of the current techniques to monitor bone condition use 
bone mass measurements. However bone mass measurements are not reliable to distinguish 
Osteoporotic and Normal subjects [29]. It has been shown that the overlap between normals 
and osteoporosis is found for all of the bone mass measurement technologies: single and dual 
photon absorptiometry, quantitative computed tomography and direct measurement of bone 
area/ volume on biopsy as well as radiogrammetry. A similar discordance is noted in the fact 
that it has not been regularly possible to find the expected correlation between severity of 
osteoporosis and degree of bone loss [30]. 

Structural parameters such as trabecular connectivity have been proposed as features for 
assessing bone conditions [29]. It has been shown that in vertebral crush fracture patients, 
elements such as vertical trabeculae are retained more or less intact, while elements such as 
horizontal bracing trabeculae are resorbed entirely [31] [32] . This results in disconnection of 
large number of trabecular elements. However, in non-fracture patients connections between 
elements were preserved. Long vertical trabeculae are subject to buckling under loading. 
When they lose their lateral connections to adjacent trabeculae, the degree of buckling may 
exceed the inherent strength of the bone. Structure can be thus be seen as an important 
feature in assessing bone condition. 


2 Significance 

Protecting humans against extreme environmental conditions requires a thorough under- 
standing of the pathophysiological changes resulting from the exposure to those extreme 
conditions. The knowledge of the degree of medical risk associated with the exposure is of 
paramount importance in the design of effective prophylactic and therapeutic measures for 
space exploration. Major health hazards due to musculoskeletal systems include the signs 
and symptoms of hypercalciuria, lengthy recovery of lost bone tissue after flight, possibility 
of irreversible trabecular bone loss, the possible effect of calcification in the soft tissues, and 
the possible increase in fracture potential. Our research to relate the local trabecular struc- 
tural information to microgravity conditions is an important initial step in understanding the 
effect of microgravity and countermeasures on bone condition and strength. The proposed 
research is also closly linked with Osteoporosis and will benefit the general population. 

3 Hypothesis/Rationale 

• The rationale for this research is based on the premise that microgravity conditions 
change bone structure as well as bone mass. 

• Bone structure can be characterized by fractal geometry. 
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• Fractal characterization of bone structural changes due to microgravity conditions is 
not only optimal but also pragmatic. 


4 Specific Tasks Achieved 

The overall goal is the characterization of bone structural changes due to microgravity with 
the aid of fractals. 

• We have related the fractal dimension computed from MRI images to those obtained 
with X-Rays of the same specimen. 

• We have computed the fractal dimension from MRI images of subjects before, during 
and after bed rest. 

• We have computed the fractal dimension of faxitron images of Osteoporotic and Normal 
Bone slices of human subjects. 

During the summer fellowship, we have implemented and tested the fractal algorithms 
on sample MRI bone images. 

The long term goal we are pursuing is to integrate fractal and finite element analysis. 
We are interested in studying the effect of various external parameters on local structural 
information provided by our method. This makes it possible to quantify the effect of coun- 
termeasures on local structural information of the bones. 


5 Summary of Results Achieved 

• Fractal Dimension of the human subject tested is higher before bedrest than during 
bedrest. This shows that the Fractal Dimension decreases due to bedrest. 

• Fractal Dimension of Osteoporotic subjects is less than the Fractal Dimension of Nor- 
mal Subjects in the Faxitron Image Study. 

• Fractal Dimension measured with the MRI image and X-Ray image of the Patella is 
the same. 

6 Fractal Model 

With natural objects, familiar metrics from classical geometry such as length, area and 
volume depend on the scale at which we decide to look at the object ( e.g . the size of the 
“yard stick”). As an example, one can show that the surface area of a sand grain is entirely 
dependent on the scale at which one chooses to look at it. The smaller the scale, higher the 
surface area since more nooks and crannies become visible. 
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Fractal geometry [12] characterizes this ability of an n dimensional object to fill the n + 1 
dimensional space where the relationship of a measure M with topological dimension n and 
the scale e is expressed as 

M(e) oc 1 0 < r < 1 (1) 

where the quantity r + n (denoted by D) is called the fractal dimension or Hausdorj- 
Besicovich dimension and characterizes the degree of erratic behavior. For an ideal object, 
the measure M is independent of the scale t and hence n = D. Thus a fractal object can be 
defined as a set for which the fractal dimension is greater than the topological dimension. 
There are number of ways in which fractal geometry can be applied to the analysis of texture 
in images. Pentland et.al [19] used a method related to the cooccurance matrix technique of 
texture classification based on fractal dimension. They find the standard deviations of the 
difference of gray levels separated by a given vector and plot it against the vector lengths as 
a log-log graph. In another technique [18] two dimensional gray level image is represented 
as a three-dimensional surface whose height at each point represents the gray level at that 
point and the surface area is measured at different scales. It has also been shown that [9] 
a n- dimensional fractal object can be characterized by the fractional Brownian motion of n 
variables and that the relationship between the power spectral density and r are independent 
of the projection [7]. This makes the fractal dimension computed from the projections of 
n-dimensional fractal object represent that of the original object. 

It has been shown that the power spectra of a fractal object exhibits an inverse-power 
relation to frequency. For one- dimensional signals, this can be expressed as, 

-S’(w) oc -p^T ( 2 ) 

M 


and the parameter r is related to the fractal dimension by the equation, 


D = 


5 — r 

2 


( 3 ) 


Thus one can exploit the above relationships to estimate the fractal dimension in frequency 
domain by finding the slope of the plot of log power spectra vs. log frequency. For multi- 
dimensional case, the corresponding relationship has been shown to be 


where 


S(u-l, UJ 2 ,- ■■ OC 


1 


( 4 ) 


_ _ 2n + 3 - r 

“ 2 

In the next section, we reiterate some basic morphological operations and present an 
algorithm to obtain the fractal dimension of a three dimensional surface. This approach gives 
us the added capability of dealing with shape using different structuring elements. Since the 
underlying theme is based on fractals, all of the above properties hold for morphological 
fractals as well. 
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Our approach is based on the above method where the digitized radiograph is represented 
as a surface whose height represents the gray level at each point. The surface area at different 
scales is estimated using a series of dilations of this surface by a given structuring element 
whose size determines the scale. The derivations in the next section shows how one can 
obtain the surface area from the volume of the dilates and some simplifications which allows 
one to use dilations by a fixed size structuring element. 

7 MORPHOLOGICAL FRACTALS 

Mathematical Morphology as developed by Matheron and Serra [25] is basically a Set Theory 
and uses set transformations for Image Analysis. It extracts the impact of a particular 
shape on images via the concept of Structuring Element (SE). The SE encodes the primitive 
shape information. In a discrete approach, the shape is described as a set of vectors with 
respect to a particular point, the Center, which does not necessarily belong to the SE. 
During Morphological transformation, the Center scans the whole image and matching shape 
information is used to define the transformation. The transformed image is thus a function 
of the SE distribution in the original image. 

In particular, Dilation of a set X with a SE Y is given by the expression 

X®Y ={x:Y*{jXj:$} (6) 

where Y x indicates the translation of set Y with x. The operations gives a set whose surface 
is the path traced by the center of the SE Y when it traverses along the surface of X. Using 
the above operation, surface area of a compact set X with respect to a compact convex SE 
Y which is symmetrical with respect to the origin is given by [25] : 


5(X,U) = lim 

v ' P-0 


V(dX ® P Y) 

2p 


( 7 ) 


where dX is the boundary of set X and © denotes the dilation of the boundary of set 
X by the SE Y scaled by a factor p. V (X) gives the volume of set X It has been observed 
that even though for “regular” classes of sets, the surface area S(X, Y) is finite, for many 
“natural” objects, this can be infinite. 

From the above expressions, it can be seen that dilating by pY hides all structures smaller 
than pY and therefore is equivalent to looking at the surface at scale p. For experimental 
purposes, we can calculate the surface area of a set X at scale p 


S(dX,Y,p ) 


V{dX © pY) 

2p 


( 8 ) 


If the object is regular, the surface area will not change with /?,. For a fractal object, 
S(dX, Y,p) increases exponentially with decreasing p as seen from equation (1). By taking 
the logarithm, we now have 
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log(S(<«,J») = log( A^) — r log(/>) 

D = 2 + r (9) 

Where K is the proportionality constant. We can now estimate D by plotting log[S(<9X, Y, p,)] 
Vs. log[p,] for a given set of scale factors pi, i = 1,2, • • • ,N and calculating the gradient of 
the line that fits the plot. 

The series of dilation of X by p t Y required for the above computation can be further 
reduced to dilation by the unit element Y by observing that if Xy = X 0 pY, then 

Xf 1 = X®(p + 1)Y = {X®pY)®Y = X p b ®Y (10) 

Apart from having projection angle and scale invariance, with the morphology method, 
we also have the freedom of selecting a structuring element suited for the problem at hand. 
This, coupled with the fact that the method only involves dilation makes the implementation 
straightforward compared with the other methods where it is required to find the covering 
of the boundary of X ( dX ) by spheres in 3D and disks in 2D. 

The surface area S(X , Y, pi) can be iteratively calculated as follows: Let the image X be 
defined as the set of triplets {< f(x , y)-,x = l,N;y = l, M} and the the structuring element 
y be given as a set of triplets {< x, , yi, Zi >, i — 1, P}. The p th dilate f p (x, y ) is calculated 
as 


f p (x,y) = max{/ p _i(r + x { ,y + yi) + z { , i - 1,2, ••• ,P] 


The initial condition f 0 (x,y) is set to f(x,y). 

The surface area at each step can be calculated by using equation(8) where 


S{X,Y,p) 


V{X,Y,p) _ Ex€l,iV;v61.M(/p(^ l y) ~/(^3/)) 

2 p 2 p 


(ii) 


(12) 


8 Alternating Sequential Filters 

In this section we will use Morphological Pyramids to compute the Fractal Dimension. Gen- 
erally, dilation and erosion are applied in pairs to make the transformation independent of 
the origin of the SE. The opening operation of a set F by a SE B is defined as, 


FoB = (FQB)®B 


The dual operation closing is defined as 


F»B = (F@B)eB 


A composite opening- closing mapping can be defined as, 


Mb{F) = (F o B) • B 
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An iterative application of such an operator is defined as a Alternating Sequential Fil- 
ter(ASF) [15], [26]. 


ASF(F) = M Btt M B „_ l ...M Bi {F) 

where N is an integer and Bn, Bn-i,.Bi are SEs with different sizes satisfying the 
constraint Bn D Bn - i 2 ~ B\. 

During a morphological transformation, the structuring element scans the whole image 
and modifies each point depending on the structural similarity of the SE with the image 
at that point. During erosion, any structure in the foreground smaller than the SE is 
removed from the image. Similarly, dilation removes such structures in the background. 
The composite operation of opening and closing provides the same result with the added 
advantage of independence from the origin of the SE. The result of such a transformation 
can be interpreted as a transformation where the details of the original image smaller than 
the SE axe removed. Thus the repeated application of such transformations using SE with 
increasing size will result in a sequence of images with decreasing details. This is equivalent 
to a representation of the image at decreasing resolution. This multiresolution representation 
is used for the surface area measurement. 

Though opening and closing are sufficient for such a purpose, use of ASF is preferred. 
ASF is more robust and introduces less distortion than individual application of opening or 
closing [15]. A smaller structuring element should be used before a larger one. 

The direct application of ASF needs iterative morphological transformation using in- 
creasing structuring elements. As the SE size increases, the computation involves increases. 
In [15], it has been shown using morphological sampling theorem [4] that ASF pyramid can 
be obtained equivalently by decreasing the size of the image instead of increasing the size of 
structuring element. At any level of the pyramid, the next level representation is obtained 
by subsampling the image and then transforming the sampled image using a constant sized 
element. As the image size is decreased and the SE size remains unchanged, this process is 
more efficient. 

Once the morphological pyramid is obtained, the surface area is computed using a piece- 
wise planer approximation. The image is viewed as a sampled version of continuous two- 
dimensional surface. In the approximation, the surface is considered as the union of non- 
overlapping triangular areas. The triangles axe formed in two steps. First, the support of the 
surface is divided into a number of squares, each having the side equal to one grid spacing. 
An M x N image is composed of MN such unit squares. Each square is represented by 
the pixels at four corners (p,q),(p+ l,q),(p,q + 1) and (p + \,q + 1). In the next step, 
each such square is divided into two triangles, one of which is represented by the pixels at 
(p, q), (p+1,9) and (p,q + 1) and the other by (p + l,$), (p+l, 9 + l) and (p,q-\- 1). The area 
of any such triangle is computed using the pixel gray level values of its three corners. The 
surface area of the image is approximated by the sum of the area of all the triangles. The 
image size decreases as the resolution decreases. To represent the correct surface area at all 
resolutions, the area computation takes the resolution step into account by normalizing the 
grid space (for example, the grid spacing in the 1/4-th resolution represents 4 gridspacing of 
the original resolution). 
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The slope of the log-log plot of the surface area against the resolution is computed using 
least-square estimation. 


9 Preliminary Experimental Results 

9.1 Relating the Fractal Dimension of X-Ray and MRI images 
of the patella. 

An isolated human patella was scanned in MRI scanner. Additionally, X-Ray images of the 
same patella was obtained. Fractal dimension was computed in four different regions of the 
MRI image with a window size of 50. Fractal Dimension was also computed in the same 
regions on the X-Ray image. Figure 1 shows the MRI image and the X-Ray image of the 
patella. The computed fractal dimension in each of the four different areas are shown. The 
tabulated results in figure 1 show that there is no significant difference between the fractal 
dimension computations between the MRI image and the X-Ray image. 

9.2 Faxitron Bone Images of Osteoporotic and Normal Subjects 

Faxitron bone images of the normal and human subjects were obtained. Fractal dimension 
was computed with a window size of 110 by 110 on six slices of normal and 8 slices of 
osteoporotic bone images. The fractal dimension computed with the flat and pyramid SEs 
are shown. The fractal dimension computed on the osteoporotic subjects was smaller than 
the fractal dimension computed with the normal subjects. 

9.3 Bedrest Studies 

Fractal dimension was computed on a series of MRI images of the vertebra of human subjects 
scanned before, during and after bedrest. The computed result on figure 3 shows that the 
fractal dimension of the subject before bedrest was higher than the fractal dimension during 
bedrest. 
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MRI image of patella 


X-Ray image of patella 



I s/patella_pieces2/MRI-patella.l. 1 50 1 D=2.799259 1 
I s/patelIa_pieces 2 /MRI-patelIa. 2 . 1 50 1 D=2.797808 1 
I s/pateIIa_pieces 2 /MRI-patel!a. 3 . 1 50 1 D=2.736876 1 
I s/patella_pieces 2 /MRI-pateUa. 4 . 1 50 1 D=2.694368 1 

I s/pateIIa_pieces2/XRAY_patella.l. 1 50 1 D=2.804438 1 
I s/pateIla_pieces 2 /XRAY__patella. 2 . 1 50 I D=2.750343 1 
I s/patelIa_pieces2/XRAY_pateIla3. 1 50 1 D=2.735890 1 
I s/patelIa_pieces2/XRAY_patella.4. 1 50 I D= 2.622544 1 


Fig 1: Fractal Dimension of X-Ray and MRI image 
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Please note that the fractal number can vary between 2 and 3 . 

The differences obtained between the groups is quite significant. 


/***** Results from 152_L4 stack using morphology *****/ 
/***** Each region-of-interest is of dim. 110x110 *****/ 


es/152 L4 1.1. 

1 

110 

Flat 

D=2. 609848 

1 

Pyramid 
D=2. 543952 

es/152 L4 2.1. 

1 

110 

D=2. 670599 

1 

D=2. 595224 

es/152 L4 3.1. 

1 

110 

D=2 .661320 

1 

D=2. 598183 

es/152 L4 4.1. 

1 

110 

D=2. 667757 

1 

D=2. 597461 

es/152 L4 5.1. 

1 

110 

D=2. 667157 

1 

D=2. 593838 

es/152 L4 6.1. 

1 

110 

D=2. 674416 

1 

D=2. 610721 


Normal 


/★**** Results from 16_L4 stack using morphology ★****/ 
/***** Each region-of-interest is of dim. 110x110 *****/ 


s/16_L4_2 . 1 . I 
s/16_L4_3.1. I 
s/16_L4_4 . 1. I 
s/16_L4_5.1. I 
s/16_L4_6.1. I 
s/16_L4_7.1. I 
s/16_L4_8 . 1. I 
s/16_L4_9. 1. I 


110 

1 

D=2. 478800 

110 

1 

D=2. 523201 

110 

1 

D=2. 518456 

110 

1 

D=2. 514385 

110 

1 

D=2. 572345 

110 

1 

D=2. 572345 

110 

1 

D=2. 450142 

110 

1 

D=2. 482301 


D=2. 417768 I 
D=2. 451674 I 
D=2. 447915 I 
D=2. 446387 I 
D=2. 481285 I 
D=2. 481302 I 
D=2. 402326 I 
D=2. 425839 H 


Osteoprotic 


Fig 2: Fractal Dimension of Normal & Osteoprotic images. 
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FD value 



L2 




10 Conclusion 


We use fractal dimension to characterize bone structure. We see that the fractal dimension 
computed with MRI images and X-Ray images of the patella are the same. We find that 
the fractal dimension of the Osteoporotic subjects to be less than the fractal dimension of 
Normal subjects. We also show that the fractal dimension of the subjects during bedrest is 
less than the fractal dimension before bedrest. 
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